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Abstract 

We present a two-dimensional (2D) generalization of the stabilized Kuramoto 
- Sivashinsky (KS) system, based on the Kadomtsev-Petviashvili (KP) equation 
including dissipation of the generic (Newell - Whitehead - Segel, NWS) type and 
gain. The system directly applies to the description of gravity-capillary waves 
on the surface of a liquid layer flowing down an inclined plane, with a surfactant 
diffusing along the layer's surface. Actually, the model is quite general, offering 
a simple way to stabilize nonlinear media combining the weakly-2D dispersion 
of the KP type with gain and NWS dissipation. Other applications are internal 
waves in multilayer fluids flowing down the inclined plane, double-front flames 
in gaseous mixtures, etc. Parallel to this weakly 2D model, we also introduce 
and study a semi-phenomenological one, whose dissipative terms are isotropic, 
rather than of the NWS type, in order to check if qualitative results are sensitive 
to the exact form of the lossy terms. The models include an additional linear 
equation of the advection-diffusion type, linearly coupled to the main KP-NWS 
equation. The extra equation provides for stability of the zero background in 
the system, thus opening a way to the existence of stable localized pulses. We 
focus on the most interesting case, when the dispersive part of the system is of 
the KP-I type, which corresponds, e.g., to the capillary waves, and makes the 
existence of completely localized 2D pulses possible. Treating the losses and 
gain as small perturbations and making use of the balance equation for the field 
momentum, we find that the equilibrium between the gain and losses may select 
two steady-state solitons from their continuous family existing in the absence of 
the dissipative terms (the latter family is found in an exact analytical form, and 
is demonstrated numerically to be stable). The selected soliton with the larger 
amplitude is expected to be stable. Direct simulations completely corroborate 
the analytical predictions, for both the physical and phenomenological models. 

PACS: 05.45.Yv, 46.15.Ff 
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1 Introduction and derivation of the model 



Localized structures, such as solitary pulses (SPs), play a dominant role in many 
conservative and dissipative nonlinear physical systems. As is commonly known, 
in conservative systems SPs are supported by balance between nonlinearity and 
dispersion [jjj, while in dissipative models, such as the Ginzburg-Landau equa- 
tions, it must be supplemented by the balance between losses and gain Q|. 

An important example of a one-dimensional (ID) model that combines con- 
servative and dissipative effects is a mixed Kuramoto - Sivashinsky (KS) - Ko- 
rteweg - de Vries (KdV) equation, which was first introduced by Benney || 
and is therefore also called the Benney equation. This equation finds various 
applications in plasma physics, hydrodynamics and other fields (J, |). SPs are, 
obviously, important objects in systems of this type however, they cannot be 
completely stable objects in the Benney equation, as the zero solution, which 
is a background on top of which SPs are to be found, is linearly unstable in 
this equation due to the presence of the linear gain (however, if the dispersion 
part of the Benney equation is large enough, the growing perturbation, moving 
at its group velocity, does not actually overlap with the SP and therefore does 
not destroy it, see Refs. M and references therein). A stabilized version of the 
Benney equation was recently proposed in Ref. Q . It is based on the KS - KdV 
equation for a real wave field u(x,t), which is linearly coupled to an additional 
linear equation of the diffusion type for an extra real field v(x,t), that provides 
for the stabilization of the zero background: 

l^t ~t~ U1M X t^xxx V x (XU XX ^f^xxxx j (1) 

v t + cv x -u x = Tv xx . (2) 

Here, a, 7 and T are coefficients accounting for the gain and loss in the it- 
subsystem and loss in the w-subsystem, respectively, and c is a group-velocity 
mismatch between the fields. 

It was shown both analytically (by means of the perturbation theory) and 
numerically in Ref. || that the system (jl|) , (|J) gives rise to a completely stable 
SP, as well as to stable bound states of SPs, in a broad parametric region. As 
a matter of fact, Eqs. (|l|) and (||) furnish the first example of a model of the 
KS type that gives rise to fully stable pulses, and they can be easily observed 
in experiment (see a detailed discussion of physical realizations of the model - 
first of all, in terms of liquid films flowing down under the action of gravity - in 
the next section). 

The liquid films and other systems in which the observation of stable SPs is 
expected (e.g., double- front flames, see below) are two-dimensional (2D) media, 
therefore a relevant issue is to introduce a physically meaningful 2D version of 
equations (|l|) and (J2j> and seek for stable 2D pulses in the generalized model. 
This is the objective of the present work. It will be demonstrated that the 
2D model which will be derived here is a generic one for a number of different 
applications. The results will directly point at a new type of 2D pulses that can 
be observed experimentally in a straightforward way. Besides that, a possibility 
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of the existence of stable 2D pulses may help to understand the phenomenon of 
turbulent spots, see, e.g., Refs. || and references therein. 

The paper is organized as follows. The 2D model is derived in detail in 
section 2, starting with a particular physical problem, viz., a downflowing liquid 
film carrying a surfactant, and then proceeding to the generic character of the 
model. Parallel to the derived model, we will also consider its counterpart which 
differs by the form of 2D dissipative terms, in order to demonstrate that basic 
results are insensitive to particular features of the model (which is relevant to 
show, even if the model is generic). In section 3, we consider stability of the zero 
solution in the 2D model, which, as well in the ID case, is a necessary condition 
for the full stability of SPs. In section 4, an analytical perturbation theory for 
SPs is developed by treating the gain and loss constants as small parameters. 
To this end, first a family of exact 2D soliton solutions of the zero-order system 
(the one without the gain and loss terms) is obtained, following the pattern of 
the well-known "lump" solitons of the KP-I equation. Then, using the balance 
equation for the net field momentum, similar to how it was done in the ID model 
[|| , we demonstrate that the combination of the gain and loss terms may select 
two (or none) stationary pulses out of the continuous soliton family existing 
in the zero-order system; in the case when two stationary pulses are found, it 
is very plausible that the one with the larger value of its amplitude is stable. 
In section 5, we present results of direct numerical simulations of the full 2D 
model, which completely confirm the analytical predictions, i.e., the existence 
of stable 2D localized SPs. The paper is concluded by section 6. 



2 The model 

The physical meaning of the model based on Eqs. ([I]) and (||), and its 2D gener- 
alization developed below, can be understood in terms of a particular application 
to a thin downflowing liquid layer with a surfactant trapped on its surface. As 
is well established (see a review ||), in the ID case the evolution of the flow 
velocity field v{x,t) in the layer is governed, in appropriately chosen units, by 
the KS-KdV (Benney) equation, it t + uu x + u xxx — —au xx — "fu xxxx . In this 
equation, the gain a is induced by gravity, and the loss parameter 7 is propor- 
tional to the fluid's viscosity coefficient, while the left-hand side is generated by 
the Euler's equations (for an irrotational flow), exactly the same way as in the 
classical derivation of the KdV equation. If the surfactant is distributed on the 
surface of the layer with a density c + v(x, t), where c and v are, respectively, its 
constant and small variable parts (\v\ <C c), the gradient of v creates, through 
the variation of the surface tension, an additional force ~ v x which drives the 
flow, hence the Euler's equation adds the term v x to the right-hand side of the 
KS-KdV equation (the coefficient in front of this term may be scaled to be 
1), so that the equation takes precisely the form ([!]). Further, the evolution of 
the surfactant density is governed by an obvious advection-diffusion equation: 
v t + [u (c + v)] = Tv xx , where T is the surface diffusion constant. With regard 
to the condition \v\ <C c, the latter equation takes the form (0). 
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In the 2D case, we consider a quasi-lD (weakly 2D) flow of the film, with 
the y-scale much larger than that along the x axis. In other words, if the wave 
is taken as exp (iKx + iQy), the small wavenumbers are ordered so that 



Q ~ K 2 . (3) 

Then, according to the classical derivation Jl0| |, the KdV part of Eq. ([!]) is 
replaced either by the KP-I equation, 

(U t + UU X + U XXX ) x = Uyy , (4) 

or by the KP-II equation, which is 

(U t + UU X + U XXX ) x = -Uyy , (5) 

the coordinate y being properly rescaled. The choice between Eqs. (f|) and (0) is 



determined by the sign of the 2D correction to the dispersion; in particular, the 
capillarity gives rise to KP-I. The difference between the KP-I and KP-II equa- 
tions is that, although both of them have quasi-lD (i.e., y-independent) soliton 
solutions that reduce to the usual KdV solitons, only in the KP-II equation this 
soliton is stable against y-dependent perturbations. On the other hand, the 
KP-I equation has stable 2D solitons (the so-called "lumps"), which are weakly 
(non-exponcntially) localized in both x and y, see below; the KP-II equation 
does not have 2D solitons. 

The next step is to accordingly generalize, for the 2D situation, the dissi- 
pative term in Eq. (|l]). Dissipative generalizations of the KP equations were 
introduced in some works, see, e.g., Refs. [[ll). Those generalizations follow the 
pattern of the arrangement of the KP equations proper: if one starts from a 
corresponding ID equation containing dissipative terms [for instance, Eq. (ffl)], 



which is written as something = 0, its 2D counterpart is (something) x = ±u yy . 
The accordingly modified Eqs. ([j]) and @ then take the form 

(u t +uu x + u xxx — v x + au xx +^u xxxx ) x = ±u yy , (6) 
(v t + cv x -u x - Tv xx ) x = ±v yy . (7) 



The second equation in this system can be simplified, as the usual ordering 
of the partial derivatives, adopted in course of the quasi-lD derivation ||lo| , 
implies that v yy is a small quantity of a higher order than v xxx , hence v yy may 
be dropped, so that Eq. ((?]) remains one-dimensional. As a result, this version 
of the 2D system takes the form 

uu x -t~ u xxx v x -\- glil xx -\- 'yu xxxx ^) x — iiAyy, (8) 
Vt + cv x -u x - Tv xx = 0. (9) 

However, while keeping Eq. (g) in the ID form is quite acceptable, the 
way the dissipative and gain terms in Eq. (||) were made two-dimensional was 
formal, being not based on any physical argument. Moreover, it will be shown 
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in the next section that, unlike its ID counterpart (|l]) and (0), this 2D model 
cannot produce any stable solitary pulses, as its zero solution is always unstable. 

In order to derive a physically relevant form of the dissipative part of the it- 
equation, one should resort to the standard procedure that derives a generic set 
of dissipative terms in the quasi-lD situation (in the context of convective flows) 
in the Newell- Whitehead-Segel (NWS) equation [Q. This equation gives a 
simple prescription, which, as well as the KP equations, is based on the ordering 
(||) of the longitudinal and transverse wavenumbers: the longitudinal dissipative 
term ~/u xxxx must be supplemented by its transverse counterpart ju yv [the 
scaling of the transverse coordinate y in the NWS equation is precisely the 
same which cast the KP equation in the standard form (||) or <Js[) ; the identity 
of the scalings is not accidental, being a consequence of the fact that both KP 
and NWS equations are generic ones in the quasi- ID geometry, the former one 
in the class of dispersive equations, and the latter among dissipative equations]. 
Thus, the proper form of the 2D system is 

(Ut + UU x +U xxx -V X +aU xx +~fU XXXX -JUyy) x = ±Uyy, (10) 

v t + cv x - u x -Tv xx = 0. (11) 

In fact, exactly the same combination of dissipative terms as in Eq. ( |l0| ) has 
been earlier derived in asymptotic equations governing nonlinear waves on thin 
downfl owing liquid films M,Jl4] and in two- fluid flows 



Note that Eqs. ( |10| ) and (11) does not contain any additional free parameter 
in comparison with the ID system (|j), (||). As a matter of fact, this is another 
consequence of the fact that both the dispersive and dissipative parts of the 
system were derived for the generic quasi- ID case. 

A special form of a 2D quasi-isotropic (rather than quasi-lD) generalization 
of the Benney equation was also derived, which, for instance, describes Rossby 
waves in a rotating atmosphere |l6| (see also Refs. [jl3[ Q): 

u t + uu x + Au x + au xx + 7 A 2 u = 0, (12) 

where A = d x + <9 2 , hence A 2 in Eq. ( p^) is a fourth-order isotropic dissipative 
operator. Following the pattern of Eq. (|12[), a quasi-isotropic generalization of 
Eqs. ([[]) and (Q) may be introduced, replacing the term u xxxx in Eq. (||) by A 2 , 
which leads to a system 

(u t + uu x + u xxx - v x + au xx + jA 2 u) x = ±u yy , (13) 
v t + cv x -u x - Tv xx = 0. (14) 

Note that all the 2D systems introduced above conserve two masses, 

/> + oo p-\-oo />+oo p-\-oo 

M= / / u(x,y)dx dy, N = / / v(x,y)dxdy (15) 



— oo J — oo 



(which indeed have the meaning of masses in the application to liquid-film flows, 
in the cases of both the single layer with surfactant and two layers with the lower 
one dominated by viscosity). 
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Thus, Eqs. @) and (JTT|) can be derived in a consistent way as a system 
describing the downflow of a liquid viscous him carrying a surfactant. Actually, 
the derivation outlined above clearly suggests that this model is a generic one 
for weakly 2D systems combining dispersion, gain, and viscosity. In particular, a 
derivation involving more technicalities and following the lines of Ref. |5| shows 
that the same model applies to a downflow of a two-layered liquid film in the 
case when the substrate layer is dominated by viscosity. A physically different 
example may be a double-front flame propagating in a combustible gaseous 
mixture, in the well-studied case when the combustion involves two consecutive 
reactions (see Refs. |l7j and references therein). In the case when a single- 
flame front is unstable, it is well known that its evolution is governed by the 
KS equation proper [Q . It is known too that a situation with one front stable 
and one unstable is possible, which is described by a linearly coupled system, 
consisting of a KS equation and the one tantamount to Eq. ([l3|). Dispersion, 
which is missing in the KS equation proper, can be induced by a background 
shear flow tangent to the flame , but detailed derivation of the full model 
for this case is beyond the scope of this work. 

As concerns the model ( |l3| ) , (|lj) , we consider it as a semi-phenomenological 
one, which may apply to cases which are "more isotropic" than those obeying 
the condition (||). We will study this model parallel to the physical one (|l0|), 
( p"l| ) in order to see if qualitative results are sensitive to details of a given model. 
Accordingly, the systems (|io|), ( p"l| ) and (13), ( p^ ) will be below referred to as 
physical and phenomenological ones, respectively. In the analysis presented in 
the following sections, we focus on the case of the KP-I type, i.e., with the upper 
sign in Eqs. (||), ( |l0| ) and (Jll|), as only in this case one may expect the existence 
of nontrivial 2D pulses, while the models of the KP-II type may only extend 
the SP found in Ref. M into a quasi-lD (y-independent) pulse in 2D. 



3 The stability of the zero solution 

As previously mentioned, completely stable SPs can only exist in a system 
whose trivial solution, u = v = 0, is stable, therefore our first objective is to 
analyze this necessary stability condition. We substitute into the corresponding 
linearized equations a 2D perturbation in the form u ~ exp (ikx + iqy + Xt), 
v ~ exp (ikx + iqy + Xt), where k and q are arbitrary real wave numbers of the 
perturbation, and A is the corresponding instability growth rate [note that, as 
all the equations that we are going to consider are written in the scaled form, k 
and q are not assumed to be specially small, unlike K and Q in Eq. (g)]. This 
leads to a linearized dispersion equation 

[k (A - ik 3 - ak 2 + jk 4 ) - iq 2 ] (A + ick + Tk 2 ) + k 3 = (16) 

for the formal model (||), (^|), or 

[k (X - ik 3 - ak 2 + 7 (fc 2 + q 2 ) 2 ) - iq 2 ] (X + ick + Tk 2 ) + k 3 = (17) 
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for the phenomenological model (|l3|), (fu]), or 

[k (A - ik 3 - ak 2 + 7 (fc 4 + q 2 )) - iq 2 ] (A + ick + Yk 2 ) + k 3 = (18) 

for the physical model (|l0|), (pi]). The stability condition states that both 
solutions of the quadratic equations ([l6]) , ( p7[ ) , or ( |l8| ) must satisfy the inequality 

Re [\{k, q)} < (19) 

at all the real values of k and q. 

As it was already mentioned, the zero solution in the formal model (||), ([)[) 
is always unstable, which can be shown as follows: in the case when k is small, 
while q is ~ 0(1), then two roots of Eq. (|l^) can be expanded as 

Ai(fc) = -ick - Yk 2 + ... , A 2 = iq 2 /k + ak 2 H . (20) 

Obviously, the second root in Eq. ( p0| ) yields instability. 

The zero solution may be stable in the physical and phenomenological mod- 



els. Although direct check of the condition (|19|) for Eqs. ( |17| ) and (18) at all 
real values of k and q with the four free parameters is a formidable algebraic 
problem, it is possible to link the stability condition for the physical system to 
that which was studied in detail for the ID system ([!]), (|2|) in Ref. Q. An alge- 
braic transformation shows that, if the condition Re [A(fc)] < holds at all real 
values of k in Eq. (8) of Ref. ||, then the inequality Re [X(k,q)] < is true at 
all real values of k and q in Eq. (|l8|), or, in other words, the stability of the zero 
solution in the ID case guarantees its stability in the 2D case for the physical 
model. However, rather than following formal algebra, it is easier to understand 
this result from Fig. 1, which shows a 3D plot of the instability growth rate 
Re A vs. k and q for a set of typical values of the parameters. Obviously, when 
k is large, the stability is secured by the higher-order dissipative term in Eq. 
( [To| ) . The most dangerous case is that when k and q are relatively small. It can 
be seen from Fig. 1 that the growth rate, considered as a function of q, attains 
its maximum value at q = if k is very small. As k increases, two additional 
local maxima of the growth rate appear at nonzero q, and they may be greater 
than the local maximum at q = 0. However, since the local maximum at q = 
becomes more negative as k increases, the two side maxima remain negative too. 
Thus, the zero solution to Eqs. ( |Toj ) and ([ll]) is stable in the same parameter 
region in which it was found to be stable in the ID system (|l|), (||) in Ref. I|. 

In the phenomenological model, the zero solution is also stable in a certain 
parametric region. However, no simple relation of the stability condition to that 
in the ID system can be found in this model. 



4 The perturbation theory for two-dimensional 
solitary pulses 

Both the physical and phenomenological models reduce to a zero-order system 
by setting a = 7 = Y = 0, while keeping an arbitrary value of c. This zero- 
order system is conservative, consisting of the KP-I equation coupled to an extra 
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linear one, 

(tit + uu x + u xxx - v x ) x = u yy , v t + cv x = u x . (21) 



Looking for a solution to Eqs. (|21j) in the form of a soliton traveling at a 
constant velocity s in the x-direction, so that 

u(x, y, t) = u(£, y), v(x, y, t) = y), with £ = x - st, (22) 

we immediately conclude that, as well as in the ID case, for such a solution we 
have 

v(t,v) = (c-a)- 1 u(t,y). (23) 

With regard to the relation (|23|), an exact solution to Eqs. (|Tj) giving a 2D 
weakly localized soliton ( "lump" ) is 

, , 24 (1 + cs - s 2 ) , , 



c — s 



where 



•f> C = \/ — : — : L w(C,z) = TT — j — —2 • (25) 



[z 2 + C 2 + 3] z 



This solution exists in the case 

1 + cs — s 2 



> . (26) 
c — a 

Note that, in contrast to the lump soliton solutions of the KP-I equation 
proper, which may only have positive velocities, the solitons ( p4| ) may also 
move in the negative direction, provided that the velocity satisfies the condition 
(p6|). In the particular case c = 0, this condition means that either s > 1, or 
-1 < s < 0. 

We have checked by direct simulations of Eqs. ( ^l|) that all the lump-type 
SP solutions are stable within the framework of the unperturbed equations ( pl|) . 
We have also simulated catch-up and head-on collisions between two lump pulses 
with different velocities. The results clearly show that the collisions are inelastic, 
hence, on the contrary to the KP-I equation, the conservative system ( ^l[ ) is 
not an exactly integrable one. An example of the head-on collision between two 
lump-type SPs with velocities s — 2.0 and s = —0.8 is shown in Fig. 2. The 
inelastic character of the interaction is obvious. 

In fact, the collision shown in Fig. 2 is not quite generic, as it is generated by 
an initial configuration in which the mass N of the v-component, see Eq. (|l5|), 
is almost zero (TV vanishes exactly if the initial velocity of the second soliton is 
—0.756 instead of —0.8, which is the case in Fig. 2). In this case, a resonance 
probably happens, resulting in generation of an additional pair of solitons, which 
is suggested by Fig. 2 [especially, by the panel 2(d)]. Far from this resonant 
typical collision (see an example in Fig. 3) is also strongly inelastic, but 
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in this case the result may be regarded, in the first approximation, as fusion of 
two solitons into one. 

To check the accuracy of the results produced by direct simulations of the 
2D equations, including those displayed in Figs. 2 and 3 and in the next section, 
the simulations were repeated increasing the number of points. For instance, 
Figs. 2 and 3 were obtained in a 2D domain of 256 x 256 points; repeating the 
same simulations with 512 x 512 points, we had obtained the same pictures, 
without any visible difference. 

Next, we restore the terms with the coefficients a, 7 and T as small pertur- 
bations, with the aim to predict, at the first order of the perturbation theory, 
which solutions will be selected from the family of the lump solitons (Q) found 
in the zero-order system (^l]). As it is suggested by the analysis of the ID sys- 
tem it is convenient to choose the relative velocity, S = c — s, as a parameter 
of the family. 

As well as in the ID case, we select stationary pulses, using the balance 
condition for the net field momentum. 



1 

P =2 



dx 



dy [u 2 (x,y) + v 2 (x,y)] , 



(27) 



which is a dynamical invariant of the conservative (zero-order) system [the 
masses ( |T5| ) cannot be used for this purpose, as they remain dynamical in- 
variants in the full dissipative models] . The dissipative and gain terms give rise 
to the following exact balance equations for the momentum: 



+ OC 



— — = I dx 



dP 
dt 



+ OO 



dy [au 2 x - j(u 2 xx + u 2 y )- Tv 2 



in the physical system (p. Of), (11), and 
ip f+00 r+oc 



dt 



dx 



dy [aul - ~/(u xx + u yy ) 2 - Tv 2 ] 



(28) 



(29) 



in the phenomenological system (|13j), dl4|) . 

Substituting the unperturbed exact solution j24|) i nto the equilibrium con- 
dition, dP/dt = 0, which follows from Eqs. ( |28| ) or (29), we cast the equilibrium 
condition in the form 



where, for both the physical and phenomenological systems, 



Ci = 



i / : / d£ 



and the constant C2 is defined differently for the two models: 



a 



(phys 



+00 

dz 

—00 j — 00 



00 



d( 



aw 



d( 2 



dw 
dz 



(30) 



(31) 



(32) 
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= I dzj dC ( ^ + ^ ) • (33) 



ac 2 " + a 

Numerical computation of the integrals yields C x = 0.34708, C^ phys) = 0.75984, 
and C , 2 Phcn ' ) = 0.85770. Then, the equilibrium equation ( |30| ) yields the following 
cubic equation for 5 (where a = a/j, T = r/7): 

S 3 + (0.465 -c)8 2 - 5 - 0.46f = 0, (34) 

S 3 + (0.415 - c) <S 2 - <5 - 0.4lf = 0; (35) 

for the physical and phenomenological models, respectively. Note that these 
equations have the same general form as the one which selects equilibrium values 
of the parameter S in the ID system investigated in Ref. ||. 

Physical roots of Eq. (p4|) or (p5|) are defined as those which not only are 



real, but also satisfy the condition (26). Physical roots select particular lump- 



type SPs, from the family of lump solitons of the zero-order system, that remain 
steady pulses in the presence of the perturbations. As well as in the ID model, 
the existence of at least two physical roots is necessary for one pulse to be sta- 
ble (the other one, with a smaller value of the amplitude, is then automatically 
unstable) Numerical solution of Eqs. and ( |35|) verifies that both equa- 
tions have exactly two physical roots in a broad parametric region, while three 
physical roots can never be found. Interestingly, the velocities corresponding to 
these two physical roots are both positive, i.e., a stable lump SP cannot move 
in the negative x-direction. 

Regions in the parameter plane (a, T) in which exactly two physical roots 
have been found are shown, for 7 — 0.05, in Figs. 4 and 5 for the physical model, 
with c = and c = — 1, respectively, and in Fig. 6 for the phenomenological 
model with c = — 1. These regions are bounded by two solid lines in Figs. 4 
through 6. Beneath the lower solid line, Eq. (Q) or ([55|) has three real roots, 
but no more than one of them is physical, and above the upper solid line, two 
physical roots bifurcate into a pair of complex ones, leaving no physical roots. 

The dashed line is the border above which the zero solution is stable in the 
ID model (|l|), (||), hence it is also stable in the 2D physical model, according 
to the results reported in the previous section. As is seen in Fig. 4, in the 
parameter plane (a, T) of the physical system with c = 0, there is no overlapping 
between regions where the zero solution is stable, and where the perturbation 
theory selects two physical solutions for the 2D pulses, but a narrow overlapping 
region is found in the same model for c = —1, see Fig. 5. It is expected that 
stable lump-type SPs exist inside this narrow region, which is confirmed by 
direct simulations, as reported in the next section. Despite the existence of 
an overlapping region in Fig. 6, we have to check the condition (^6|) for Eq. 
(pT|), i.e., the stability of the zero solution in the phenomenological model, point 
by point, since the stability of the zero solution in the ID model does not 
guarantee the stability of the zero solution in the phenomenological system, see 
the previous section. 
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The same analysis can be performed for values of the loss parameter 7 other 
than 0.05, which was fixed in the above consideration. The results show that 
the variation of 7 produces little change in terms of the expected SP stability 
region. On the other hand, the group-velocity mismatch c affects the stability 
region significantly. We had numerically found that there is no stability region 
in both physical and phenomenological model for c > 0. Then, by inspecting 
the region of negative c, it was found that there is a critical value c cr such that 
stable pulses become possible for c < c cr . In the physical model, c cr « —0.8. 

To conclude the analytical consideration, it is necessary to stress that the 
stability conditions obtained here are only necessary ones. Obtaining a full 
set of sufficient stability conditions might be possible within the framework 
of a rigorous spectral analysis of the full system linearized around the pulse 
solution, which was done (in a numerical form) for the ID Bcnney (KS - KdV) 
equation in Refs. |Q. Extending this analysis to the two-component 2D case is 
an extremely difficult problem, which is beyond the scope of the present work. 
Nevertheless, numerical results reported in the next section strongly suggest 
that the necessary stability conditions, that were obtained above by means of 
the simple analytical methods, are, most probably, sufficient for the stability of 
the 2D pulses. 

5 Direct simulations of the two-dimensional soli- 
tary pulses 

To check predictions produced by the above analysis for both the physical and 
phenomenological models, against direct numerical simulations, we have up- 
graded the implicit pseudo-spectral scheme, which was developed previously for 
the ID model (|l|), (||) (see Appendix in Ref. M), to a scheme capable of dealing 
with both 2D models. It it relevant to mention that results obtained with a 
rather coarse grid, Ax = Ay = 0.5, and relatively large time steps, At = 0.1, 
remained virtually unchanged if reproduced with an essentially finer grid and 
smaller time step. The initial conditions were taken as the lump-soliton solutions 
of the zero-order system, but with arbitrarily modified values of the amplitude, 
in order to check whether strongly perturbed pulses relax to stable ones, i.e., 
whether the stable pulses are attractors. 

Collecting data produced by the systematic simulations, it has been found 
that, for the physical model, stable lump- type SPs do exist and are stable every- 
where inside the stability region in the (a, T) parameter plane which was pre- 
dicted by the analytical perturbation theory, see Fig. 5. Moreover, all the stable 
pulses were found to be strong attractors indeed. A typical lump-type SP with 
a = 0.2, 7 = 0.05, c = -1.0, and T = 0.55 is displayed in Fig. 7. The initial- 
state pulses selected as mentioned above definitely relax to this single stationary 
lump-type SP, provided that the initial amplitude Aq of their u-component ex- 
ceeds 1.5. For instance, starting with the initial amplitudes A = 1.96 and 
A = 13.33 [the corresponding values of the unperturbed-solitons' velocity are 
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s = 0.8 and s = 2.0, see Eq. (p5)], a lump-type SP develops with the values 
of the amplitude A num = 7.11 and A num — 7.16, and velocities s num = 1.32 
and s num = 1.33, respectively, by the moment of time t = 400. Meanwhile, the 
analytical prediction for the amplitude and velocity of the presumably stable 
pulse, given by Eq. for the same values of the parameters, is A ana \ w 7.38 
and Sanai ~ 1-35. Thus, the numerical results match the theoretical prediction 
well. 

On the other hand, if the initial amplitude is too small, e.g., Aq = 0.89 
(s = 0.7), the pulse decays to zero, which is natural too, as the stable zero 
solution has its own attraction basin. Note that for the second (smaller) steady- 
state pulse, which is expected to play the role of the separatrix between the 
attraction basins of the stable pulse and zero solution, the perturbation theory 
predicts, in the same case, the amplitude A ana \ ~ 1.26, so it seems quite natural 
that the initial pulses with Aq = 1.96 and Aq = 0.89 relax, respectively, to the 
stable pulse and to zero. 

For the phenomenological model, no definite stability region for the pulses 
has been found above; instead, the stability of the zero solution at each point 
inside the region bounded by the solid lines in Fig. 6 should be checked sep- 
arately. For instance, the zero solution is found to be stable at the parameter 
point (a = 0.2, T = 0.55) in Fig. 6. At this point, the two physical roots of Eq. 
( |35| ) are 2.073 and 1.786. The amplitudes of the corresponding solitons predicted 
by the perturbation theory are A ana i w 4.72 and A anal w 1.81, respectively; the 
one with the larger amplitude may be stable according to the general principle 
formulated above. Simulations demonstrate that, indeed, there is a stable 2D 
pulse with a shape very close to the predicted one, all the initial states in the 
form of the conservative-model lump, whose amplitude exceeds 2.0, relaxing to 
this stable SP. Figure 8 shows the shapes of the u- and w-components of the 
SP at t = 400, generated by the initial configuration taken as the soliton of 
the zero-order system (^lj) with the amplitude 13.33. The amplitude of the 
thus obtained u-components is about 5.14, the discrepancy with the analytical 
prediction being less than 9%. 

It was also found that all the initial lump pulses whose amplitude is less than 
1.6 decay to zero. This fact complies well with the expectation that the second 
(smaller-amplitude) steady-state pulse, predicted by the perturbation theory, 
whose amplitude is 1.81, ought to be the separatrix between the attraction 
basins of the stable pulse and zero solution in the phenomenological system. 

The simulations also show that, even if the zero background is unstable, 
lump-type SPs may be quite stable when the integration domain (supplemented 
by periodic boundary conditions in both x and y) is not very large. An explana- 
tion given in Ref. || to a similar "overstability" effect found in the ID model, 
which is based on suppression of nascent perturbations by the pulse periodically 
traveling around the domain, applies to the 2D case as well. 

Numerical simulations of interactions between two or more stable lump-type 
SPs were also carried out. Unlike the ID system (|l|), (||), in which stable bound 
states of two and three pulses were easily found M (and unlike a physically 
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relevant 2D model of a different type, which follows the pattern of the Zakharov- 
Kuznetsov equation ]2(J and will be considered elsewhere), bound states of 
pulses have not been found in both the physical and phenomenological systems 
based on the KP-I equation. Actually, this fact may be benign for experimental 
observation of the 2D pulses, as there will be no competition with multi-humped 
structures that might render the picture much more complex. 

6 Conclusion 

In this work, we have extended the ID stabilized Kuramoto - Sivashinsky sys- 
tem to the 2D case. The 2D model is, quite naturally, based on the Kadomtsev- 
Petviashvili (KP) equation supplemented by loss and gain terms. A model with 
the dissipative terms of the Newell- Whitehead-Segel type was derived in a con- 
sistent way for the weakly 2D case. The derivation was outlined for a particular 
physical system, a downflowing liquid film carrying a surfactant that diffuses 
along its surface. It was argued that, in fact, the derived model is generic, 
therefore it also applies, for instance, to the description of interfacial waves in 
a two-layered flowing film, and double-front flames in a gaseous mixture. An- 
other model, which may apply to a more isotropic situation, was put forward 
as a semi-phenomenological one. Both models consist of the generalized KP- 
I equation with gain and loss terms, linearly coupled to an extra equation of 
the advection-diffusion type. The additional linear equation stabilizes the sys- 
tem's zero solution, thus paving the way to the existence of completely stable 
2D localized solitary pulses, which are objects of an obvious physical interest. 
Treating the losses and gain as small perturbations, and employing the balance 
equation for the net field momentum, we have found that the condition of the 
equilibrium between the losses and gain may select two steady-state 2D solitons 
from their continuous family, which was found, in an exact analytical form, in 
the absence of the loss and gain (the exact solitons are similar to the "lump" 
solutions of the KP-I equation). When the zero solution is stable and, simul- 
taneously, two lump-type pulses are picked up by the balance equation for the 
momentum, the pulse with the larger value of the amplitude is expected to be 
stable in the infinitely long system, while the other pulse must be unstable, 
playing the role of a separatrix between attraction domains of the stable pulse 
and zero solution. These predictions have been completely confirmed by direct 
simulations for both the physical and phenomenological systems. Another note- 
worthy finding is that, unlike their ID counterpart, both KP-based 2D systems 
do not generate stable bound states of pulses. 
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Figure Captions 



Fig. 1. The instability growth rate Re A for the zero solution in the physical 
model [Eqs. ( |To| ) and ([H])] vs. the longitudinal and transverse perturbation 
wave numbers k and q. The values of the parameters are a = 0.2, 7 = 0.05, 
c = -1.0, and T = 0.55. 

Fig. 2. An example of the inelastic collision between two stable lump-type 
solitons with different velocities, simulated within the framework of the zero- 
order conservative system ( f2l| ) with c = 0. The u- and u-fields prior to the 
collision (at t = 0) are shown in panels (a) and (b), and after the collision (at 
t = 40) - in panels (c) and (d). The initial velocities of the two solitons are 
S! = 2.0 and s 2 = -0.8. 

Fig. 3. The inelastic collision in the conservative system (^l|) with c = 
in the case when the initial velocities of the two solitons are Si = 2.0 and 
S2 = —0.5. The panels have the same meaning as in Fig. 2, with a difference 
that the final state is displayed for t = 45. 

Fig. 4. The stability region of the zero solution, and the region inside which 
two physical roots exist for Eq. ([m]) with 7 = 0.05 and c = 0. The zero solution 
is stable above the dashed line, and inside the region bounded by two solid lines, 
Eq. ([li]) produces two physical solutions. 

Fig. 5. The expected stability region for the solitary pulses in the parametric 
plane (a,T) for 7 = 0.05 and c = —1.0. The solid and dashed lines have the 
same meaning as in Fig. 4. 

Fig. 6. The stability region of the zero solution of the ID model (|l|), (||) and 
the region inside which two physical roots exist for Eq. ( |35| ) when 7 = 0.05 and 
c = — 1 . The solid and dashed lines have the same meaning as in Fig. 4. 

Fig. 7. A stable lump-type solitary pulse in the physical model, found as 
a result of long numerically simulated evolution in the case a = 0.2, 7 = 0.05, 
c = —1.0, and T = 0.55. The panels (a) and (b) show established shapes of the 
u and v fields, respectively. 

Fig. 8. A stable lump-type solitary pulse in the phenomenological model, 
found the same way as the pulse shown in Fig. 7 for the same values of the 
parameters. The panels (a) and (b) have the same meaning as in Fig. 7. 
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